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Abstract 

We derive a series of moment equations describing the motion and flavor trans- 
formation of neutrinos in supernova. We find a particular series of moments of 
neutrino density matrix in supernova. The emission angle distribution of neutri- 
nos is described by this series of moments. We expand the equation of neutrinos 
using these moments and obtain moment equations. We find that these moments 
have very good property of convergence and the infinite series of equations can be 
truncated to equations with a small set of moments. Using a small set of moment 
equations the required computational power is reduced by about two orders of mag- 
nitude compared to that in multi-angle simulation. The study on non-linear flavor 
transformation of neutrinos is substantially simplified using these equations. Two 
flavor system of neutrinos is also considered and new equations describing the flavor 
polarization vectors of neutrinos are found. 

PACS: 14.60.Pq, 97.60Bw 

1 Introduction 

It is realized that neutrino flavor transformation in supernova is a very complicated prob- 
lem. In addition to the vacuum oscillation effect and the effect of neutrino refraction with 
ordinary matter (Mikheyev-Smirnov-Wolfenstein effect), it is well known that neutrino- 
neutrino refraction can also be important in supernova. Neutrino self-interaction can 
be important in flavor oscillation when neutrino density is sufficiently large. In core 
collapse supernova neutrino density can be so large above the neutrino sphere that the 
non-linear flavor oscillation caused by neutrino-neutrino refraction can dominate neutrino 
flavor transformation. 

Researches on non-linear neutrino oscillation in supernova have been done by many 
groups. Many interesting phenomena caused by neutrino self-interaction, such as syn- 
chronized oscillation, bipolar oscillation and spectral split, have been found. Detailed 
numerical analysis and qualitative analysis have been done to understand these phenom- 
ena. An incomplete list of researches on this subject is [Il[2l[3llll[5l[6|,[71[8l[9], HOl [TT[ 
[T2l[T3l[Tl[l5l[T6lll71|18l[T9l|20l[a^ Ahhough the development in this held is very fast 
there are still many problems not answered and the present situation of research is not 
satisfactory. 
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One problem hard to explore is how neutrino flavor transformation depends on the 
emission angle distribution of neutrinos on the neutrino sphere. It is known that if all 
neutrinos are emitted in radial direction neutrino self-interaction vanishes and non-linear 
flavor transformation disappears. Angular distribution of neutrino emission is essential in 
non-linear flavor transformation and a complete analysis of neutrino flavor transformation 
has to take it into account. However in the present formulation of the problem it is hard 
to study the effect of emission angle distribution of neutrinos. Previous researches on this 
problem use many angle bins in numerical simulation. It is very complicated and does 
not give us insight to physics in it. 

In the multi-angle simulation flavor evolutions of neutrinos on trajectories of all angle 
bins should be studied [1]. It is found that if Li, = 10^^ erg/s more than 500 angle bins are 
required in order to make the simulation converge. Even more angle bins are required for 
larger neutrino luminosity. Note that the energy range is also made discrete in numerical 
works. It is found that a complete numerical computation requires following the evolution 
governed by more than a million equations. It is terribly complicated and makes us hard 
to understand the physics in it. Another problem in multi-angle simulation is that exact 
spherical symmetry is essential in the simulation. It is unable to imagine how to simulate 
neutrino flavor transformation without the spherical symmetry in this approach. A better 
formulation of the problem is required to make it easier to solve numerically and to be 
understood qualitatively. 

The present scheme to study neutrino evolution is not only too complicated in practice 
but also conceptually not necessary. It is not necessary to know the detailed evolution 
of wave function of neutrino in every trajectory. We just need the flux, flavor content, 
angular distribution and energy distribution of neutrinos at given time and given position 
in space. It is enough using these information to study the flavor conversion rates of 
neutrinos as functions of radius r. It is also sufficient using these information to compute 
the energy deposit of neutrinos to plasma environment. Effect of the emission angle 
distribution of neutrinos can be carefully taken into account in moment expansion of 
neutrino distribution in phase space. The whole set of equations can be tremendously 
reduced if appropriate set of moments are found. 

In this article we present a new formulation describing the transport and flavor evolu- 
tion of neutrinos above the neutrino sphere in core-collapse supernova. Instead of consid- 
ering the neutrino flavor evolution on trajectories we study the density matrix of neutrino 
at given time and position. We do a moment expansion of the Liouville equation of 
neutrino using a series of moments of density matrix. The problem of neutrino flavor 
transformation in supernova is tremendously simplified by observing that this series of 
moments have very good property of convergence. Numerical computation on neutrino 
flavor transformation is greatly simplified using this series of moment equations. This 
series of moment equations is obtained when we examine the problem using spherical 
coordinate and find a series of functions which have very good convergence property. It 
is shown that the emission angle distribution of neutrinos can be systematically studied. 

In section |5] we rewrite the Liouville equation using spherical coordinate. In section [3] 
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Figure 1: Geometric picture of angles of the neutrino momentum intersect with f. 



we simplify the equation of motion with the assumption of spherical symmetry. In section 
m we introduce a series of moments of neutrino density matrix and motivate the use of 
them by examining their property of convergence. In section [5] we expand the Liouville 
equation of neutrino and work out a series of moment equations. In section |6] we discuss 
stationary approximation and further simplify the equations. In section [7] we discuss the 
truncation of the moment equations. In section [H] we present some numerical analysis. 
In section M "we work in two flavor system and derive equations governing the pendulum 
motion of neutrinos in flavor space. We conclude in section [TOl 

2 Equation in spherical coordinate 

In this section we derive the equation of motion of neutrinos in spherical coordinate. The 
neutron star core of the core-collapse supernova is located at the origin of the coordinate 
system. 

Consider the density matrix Pf{t, x) of neutrino with a given momentum p at position 
X and at time t. pf{t, is a 3 x 3 matrix when considering three flavors of neutrinos and 
is a 2 X 2 matrix when considering two flavors of neutrinos. For anti-neutrino we express 
the density matrix as p^. The evolution of neutrino flavor is described by the Liouville 
equation [22|[23]: 



and Hq is the Hamiltonian for vacuum oscillation, L = diag{nf., n^, ra,-} in the flavor base 
is the matter term given by charged lepton number densities rie^n^r- Gp is the Fermi 



-0 = -t[Ho + V2GfL + V2GfD, Pp.] 



(1) 



where D describes neutrino self-interaction and is expressed as 




(2) 
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constant. The right-hand side of Eq. ([T]) is a commutator of the effective Hamiltonian 
and the density matrix. For anti-neutrino the equation is similar except replacing Hq by 
—Hq-.Hq —Hq. 

The right-handed side of Eq. ([T]) governs the quantum evolution of neutrino flavor. It 
does not include the effect of streaming of relativistic neutrinos. Effect of the motion of 
neutrinos is implicitly given in the left-hand side of Eq. ([1]) which can be written as 

dt ~ dt dt dx^ dt dp^' ^ ' 

The second term in Eq. ([3]) is caused by the streaming of neutrino in space and the 
third term is caused by the time dependence of neutrino momentum. For free streaming 
neutrino the third term is zero in Cartesian coordinate. It does not vanish in spherical 
coordinate: x = {r,6,ip). Fig. [1] demonstrates that 6p, the angle of neutrino direction 
intersecting with the radial direction, changes as neutrino propagates. To figure out the 
left-hand side of Eq. ([1]) we need to know ^ and ^ . 

The spherical coordinate system is shown in Fig. [2l At position x we introduce a local 
coordinate system spanned by three orthogonal unit vectors: f, 9 and (p: 

r = sin 6 cos cpx + sin 6 sin cpy + cos 6z, (4) 

9 = cos 9 cos (px + cos 9 simpy — sin 92, (5) 

(p = — sinpx + cos py, (6) 

where x, y and z are three orthogonal unit vectors in Cartesian coordinate. They are 
shown in Fig. [21 The momentum p of a neutrino at position x can be projected to f, 9 
and p directions. 9p and pp are defined as: 

p-r = cos 9p, p ■ 9 = sm9pCospp, p ■ p = sm9psmpp, (7) 

where p = p/\p\ is a unit vector. In Appendix we compute the projection of momentum 
to three directions. Using Eq. ( l86l) we get 



dr p'' 

— = — = COSb'- 

dt 



— = — = —— sm 9p sm pp, (10) 
dt p" r sm 9 



where p^ is the zero component of momentum. 
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Figure 2: Spherical coordinate and the local coordinate at x = {r,6,ip). 
^ is given by the geodesic equation: 

% = -r}. P^p'/p', (11) 

where F*^ is the Christoffel symbol and summation over repeated indices is assumed. 
Equations of p"^, p^ and p"^ are computed in the Appendix. It is shown that \p\, the 
magnitude of neutrino momentum, is a constant. This is an expected result when gravi- 
tational potential is neglected. Hence in Eq. ([3]) we just need to include effect of angular 
dependence of neutrinos. Writing p0(t,x,p) = pf(t,x, \p\, cos 9p,ipp) (similarly for p^) we 
just need to include d cos 9p/dt and dipp/dt in Eq. Q. They are given in Eqs. fl57|) and 
flHSjl in Appendix. To summarize we arrive at 

— = — + cos Op—- + - sm Op cos ^p^ H —. sm dp sm <~Pp^ 

dt ot or r aO rsma dip 

+ -(1 - cos^6'p)- — ^- ctg6 sin Opsmipp—^. (12) 

r o cos Op r oipp 

3 Equation with spherical symmetry 

In the remaining part of the present article we assume that the system in core-collapse 
supernova is spherically symmetric and the angular distribution of neutrino momentum 
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at position x = (r, 9, ip) is cylindrically symmetric, that is 

dL dL 

We obtain 

+ ^os9,^-i + 1(1 - cos^^,)^ = -^[H, + V2G,{L + D),p,]. (15) 



Alternatively we introduce 



Eq. (ITS!) becomes 



Hp = 1 — cos 9p. (16) 



^+{1- f^p)^ - ^/ip(2 - = -^[^0 + V2Gp{L + D),p^]. (17) 

Pp = Ppity^: \p\:Pp) is the density matrix at time t and radius r. \p\ and Pp are variables 
used for describing distribution of neutrinos in momentum space. Pp takes value in the 
range [0, 2]. 

We note that for a realistic distribution of supernova neutrinos the density matrix 
Pp- should be a smooth function of |p| and pp. Moreover, should decrease to zero as 
Pp —>■ 1, i.e. as 9p —>■ 7i/2. This is because no neutrinos inside the neutrino sphere can be 
emitted in the direction with 9pQ > 7r/2 on the neutrino sphere where 9pQ is the emission 
angle on the neutrino sphere as shown in Fig. [H Q Neutrinos moving out of the supernova 
have < 7r/2 for r > tq. This is expressed by the following condition for the density 
matrix 

pjj{t, r, \p\,Pp) = 0, for /ip > 1 and r > Tq. (18) 

This condition will be used when we do moment expansion of equation. 
A simple application of Eq. ffTSjl is for the total intensity of neutrino 

fp = Tr[p^], (19) 

where Tr is the trace of matrix. The right-hand side of Eq. f[T^ is a commutator of 
matrices and its trace is zero. We find 

^ + cos9p^ + -(1 - cos^ Op)-^ = 0. (20) 
ot or r o cos 9p 



*0pO = "■72 at r = To means that it is a trajectory passing through the surface of the sphere. No 
neutrino can be emitted in this direction on the neutrino sphere. 
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We integrate Eq. (l20l) by J dQj^ which is over all solid angle. Since the integration limit 
does not depend on r, we can change the order of ^ and the integration symbol J dVt^. 
After integration by part using condition Eq. (fTSll we find 



dt dr r ' 



where 



fl = j dQp^ fp^, = j dQp^ cosOp U. (22) 

The simplest case is that neutrinos are all emitted in radial direction at neutrino 
sphere. Hence (3 = 1 and Eq. fl2Tl) is simplified as 

^ + -^ir^ = o- (23) 

at r'^ or 

The solution of this equation is 

.2 

_0 ^0 



/.°(^'0 = %°(t-(r-ro),ro). (24) 



tq is the radius of neutrino sphere. The meaning of this solution is obvious. Intensity of 
neutrino at time t and radius r is directly related to the neutrino intensity on neutrino 

2 

sphere at time t — (r — ro). Intensity is suppressed by the geometric factor ^, as expected 
in a spherically symmetric system. 

An interesting observation based on Eq. (12^ is that if emission of neutrino on neutrino 
sphere (at r = tq) has very weak dependence on time, ^ = is obtained for all r > rg. It 
means Eq. (1241) gives a stationary distribution of neutrino density in space. We denote the 
approximation ^ = as the stationary approximation. It is very useful for simplifying 
the moment equations. We come back to this approximation later. 



4 Moments of density matrix 

In this section we examine in detail the geometric characteristics of neutrino transport 
in supernova. We introduce moments of density matrix and check their convergence 
properties. We motivate moment expansion of Eq. ( ITTl) . 

Consider a neutrino emitted at neutrino sphere at position xq = (ro, ^o, V'o) with an 
emission angle 6pQ intersecting with the radial direction. It is shown in Fig. [1] At time 
t neutrino arrives at position x = (r, 6, ip) with an angle 6p intersecting with the radial 
direction f. From Fig. [T]one can see that 

sin 9p = — sin 9po. (25) 

^More discussion on this point is given in section [5J 
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From this we get 



cos^p = yi--|sin2^po (26) 

Neutrinos are more peaked in radial direction if r is larger. We can see that at radius r 
the condition Eq. (ITS!) says that 



pjy = 0, for sin^^p > — or cos6'p ^ 1 ^- (27) 

Alternatively Eq. (1271) is re-expressed as 

Pp- = 0, for fip > S{r) at radius r, (28) 
where S{r) is the geometric scaling factor: 

2/2 

S{r) = 1-^1 -rl/r^= "'rf-^ ,- (29) 
" 1 + a/1 - r^/H 

S{r) scales as r^^. 

We introduce fcth moment of the density matrix: 



p,(t,r,p = IpI) = j dnpfi';ppit,r), k = 0,l,2---. (30) 

The integration is over all solid angle of neutrinos. It's easy to see that pk converge very 
fast as neutrinos propagate out of supernova. Note that J_^d cos 9p = dpp and using 
(|28|) we find 

1-2 i'2-K fS{r) p2-K 

Pk{t,r,p)= dpp dLppPpPf;= dpp d'^p pI Pp. (31) 

JQ Jo Jo Jo 

Because Pp = for pp > S{r) the integration over pp is effectively restricted to the range 
[0, S]. Pp factor in the integration gives a strong suppression. The larger the number k 
is, the stronger the suppression is. From Eq. fl^ we note that pk is suppressed roughly 
by S''+^: 

Pk oc S''^\ (32) 

Po scales as S{r), i.e. approximately r~^, if neglecting effect of flavor conversion. This 
is the geometric factor shown in Eq. fl24|) for the total intensity. Comparing with po 
Pk is further suppressed by the factor p^ in the integration and scales as S'^^^i^r), i.e. 
approximately r~'^^''^^\ Hence at large r higher moments can be safely neglected. It's 
clear that this set of moments has very good convergence property. If we use this set of 
moments to expand Eq. ( JT71) we should have a series of equations which has very good 
convergence property. 

Moments given in Eq. fl30l) are very nice quantities supernova offered to us in studying 
flavor conversion of neutrinos. In section [5] we will use this series of moments to expand 
Eq. (HZD. 
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5 Moment expansion of the equation of motion 



In this section we derive moment equations describing the transport and flavor transfor- 
mation of neutrinos in supernova. 

First, we examine the right-hand side of Eq. (fT7|l and express the self- interaction 
term, in terms of pk- Notice that Eq. (fl^ says 



j d^q (1 - cosp-q-) = j p^ (1 - cos cos6'p). 



(33) 



where 6q and Op are angles of q and p separately. Using 

1 — cos 6*^ cos 6'p = (1 — cos6'q) + (1 — cos^^p) — (1 — cos^^q)(l — cos^^p), (34) 
we can get 

D = (no - rao)/ip + (rii - ni)(l - /ip), (35) 



where 



no(t,r) = / 7^ P,— I^T Po(t,r,|g1). (36) 



(27r)3 ; (27r)3 



niit^r) = y (1 -cos ^9)^9 = y 1^ Piit.r,\^). (37) 

Similarly for n and p. 

We integrate Eq. f|T7j) using J dVL^ /^p = /g c^V'p /q dpp p^ which is over all solid angle 
of neutrinos. Since the integration limit does not depend on r we find that j dVl^; p^ §:;^Pp — 
^(/ dVtp Pp pp-). We do integration by part using the condition Eq. (fT8|) which is 
universal for all r > tq. We get 

dg ^ dip, -pu^r) ^ + _ + 2)p,+i] = -i[Ha^ p,] - i[Hb, p^+i], (38) 

where 



and 



Ha = Ho + V2Gf{L + D^), Hb = V2GfDo, (39) 



^0 = ^0 - ^0 - ('^i - ^1), Di = ni-ni. (40) 



Since vanishes in some region effectively there is an integration hmit which depends on r as 
shown in Eq. (PT|) . This effective Hmit does not change the fact that the differential and integral 
operators can be interchanged. It can be checked as follows. ^(/ dp,p fi^ pjj) = ■§p{J(f''^^ dfip p^ pp) = 

!jfp l4 §:Pp + (^^Mp Pp-)Up=s(r) = Jdp.p /ip S^Pp- In the last step Ppir, fip = S{r)) = from Eq. 
has been used. 
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In the remaining part of the article we will also use 



Hc = Ha + Hb (41) 

Eq. ( l38ll can be rewritten as 

IT + - = -i|ff.,P.l -i|ffB,P.«l. (42) 

For /c = and k = 1 Eq. fj42l) gives 



(9po , 1 9(rVo) 1 d{r^pi) 

-Q^ + - -i—Q^ = -AHA,P,]-^[HB,Pli (43) 

1 9(rVi) 1 9(rV2) .rr. 1 -rrr 1 /..^ 

^ + 7^^^ - = -^[HA,P^]-^[Hs,P.]. (44) 

If neglecting other terms the first two terms in left-hand side of Eq. (1431) make po scales 
as which is the correct geometric factor as shown in Eq. fl32|) . The presence of pi says 
that the average velocity of neutrino gas is smaller than the speed of light. It modifies 
the scaling behavior of po- The first two terms in the left-hand side of Eq. (jH]) give the 
geometric scaling factor as observed in Eq. ( I32i) . The presence of p2 in the left-hand 
side of Eqs. (jH]) modifies the scaling behavior of pi. Similarly, the left-hand side of Eq. 
(!38|) makes pk scale as r~'^^''~^^^ if neglecting p^+i and the presence of pk+i modifies the 
scaling behavior. 

Eq. (l38l) systematically takes into account the effect of the angular distribution of 
neutrino emission. The effect is encoded in moments introduced in the present article. In 
particular, the scaling law of the strength of the effective Hamiltonian is modified when 
including higher moments pk with k > 2. Furthermore, this series of moment equations 
has very good property of convergence and is better for studying flavor transformation of 
neutrino in supernova, for studying non-linear flavor transformation in particular. This 
virtue will help us to truncate the infinite series moment equations and will be further 
discussed in section [71 

An important point when studying effect of self-interaction in flavor conversion is that 
the minimal set of equations should include equations for po and pi, Eqs. (143|) and (H^ . 
because if neutrinos are all emitted in radial direction, that is p^ = for A; > 0, neutrino 
self-interaction vanishes. 



6 Stationary approximation 



In solving the equations of neutrino evolution in supernova we assume the stationary 
approximation: 
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This approximation is assumed by noticing that the neutrino self-interaction is important 
in the region of radius less than hundreds kilometers in which neutrino density is suffi- 
ciently large. Neutrino travels through this region in ~ 10~'^s. This time scale is much 
shorter than the time scale of neutrino emission in supernova which is around 10s. There- 
fore, we can assume that neutrino luminosity and energy do not change in the time scale 
as short as around ~ 10~^s. Furthermore, we can also assume that the distribution of 
ordinary matter in supernova does not change in such a short time scale. This assumption 
is valid if the disturbance of ordinary matter propagates with a speed much smaller than 
the speed of light. 

An example supporting this assumption is given by examining neutrino intensity given 
in Eq. fl2^ . According to the formula it is easy to see that if ^ = is assumed at the 
neutrino sphere where r = tq, ^ = is satisfied at r > tq. 

We introduce 

„' _ ^2{fc+l) -/ _ ^2(fc+l) - (Af.\ 

Pk - ^ Pk, Pk - ^ Pk I4bj 

where 

r 

z = -. (47) 
ro 

p'^ and correctly take into account the geometric scaling factor. 

With the assumption of stationary approximation for the density matrix we get equa- 
tion 



Using Eq. (HSi) repeatedly we can get 



(48) 



n— 1 n—1 

-1 

1=1 1=1 



where n > 2. For n = 1 the second line in Eq. fHOl) disappears and the equation reduces 
to Eq. (148|) . Equation of p' can be similarly obtained for anti-neutrinos. 



7 Truncation of moment equations 

In practical computation we have to truncate the infinite series of moment equations. If 
magnitude of higher moments are much smaller than Pq i we can assume that p^ = 
for k > N where is an integral. In this approximation we get a set of 2(A^ -|- 1) 
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equations for pk{pk) or p'l^ip'k) where k = 0,1 ■ ■ ■ , N. We denote this approximation as 
P/v approximation. In P/v approximation we get 



^ =-r,'Ql-^[HA,p'k]-^[Hc,Ql], (50) 



where /c = 0, 1, ■ ■ ■ ,N and 



N N 



Z=fc+1 «=fc+l 

Q]^^ = 0. Similarly we have equations for pj,. 

An appropriate choice of truncation depends both on the initial condition and on the 
property of convergence of these moments. It is interesting to note that higher moments 
are naturally suppressed at the neutrino sphere. For example, 

Pkit^ro) = -j-^PQ{t,rQ) (52) 

is obtained using Eq. (!30|) if neutrino is uniformly emitted with respect to the emission 
angle OpQ. [§ pio is about 10 times smaller than po- For a practical distribution which 
may have sharp falloff at an angle close to 9po = tt/2{ppq = 1), higher moments are 
further suppressed. As an example, suppose the distribution is uniform in the region 
< Ppo < 0.9 and drops to zero rapidly at Ppo = 0.9 we find that 

Pio{t,ro) ^^0.9^Vo(t,ro)^3%po(t,ro), (53) 

Pio{t,ro) ^ A0.99 pi(t,ro) ^7%pi(t,ro) (54) 

pio has magnitude smaller than po and pi, the major quantities in the evolution problem. 

Other examples of distribution can also be checked. Considering distribution pro- 
portional to 1 — Ppo or Ppo(l ~" Ppo) which are zero at = 1 we find for these two 
examples 

P.(t,ro) = ^^^^^^Po(t,ro), (55) 
P.(t,ro) = ^^^^^^Pi(t,ro). (56) 



§As noted previously, for realistic distribution of neutrinos Pp{ro) should decrease to zero as 9pQ — > 7r/2. 
Uniform distribution can be taken as an approximation to a distribution which is uniform in the range 
< /ipo < 1 — e and decreases to zero rapidly in the range 1 — e < /j,po < 1 where e ^ 1 is taken as a 
small positive number. This is a good approximation for not too large k. 
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or 



Pk{t,ro) 



6 



(57) 



{k + 2){k + 3) 



12 



(58) 



{k + 2){k + 3) 



We find that moments of order larger than 10 are suppressed in these models of neutrino 
emission. The precise number of moments required in study depends on the model of 
emission angle distribution of neutrinos on neutrino sphere. 

The idea that the series of moment equations can be truncated to a small set of 
equations is further supported by examining the geometric scaling property of higher 
moments. For example, according to geometric considerations p2 and are suppressed 
by factors about 1/3^ and about 1/3"^ compared to pi at r = Stq (~ 30km) if assuming 
Pk ~ Po on the neutrino sphere. A numerical study presented in section [H] shows that P5 
approximation is already quite good in the model described by Eq. fl55|) . More numerical 
analysis on the truncation of moment equations will be presented in other publications. 
We expect that a Piq approximation is probably enough to describe many phenomena in 
the transport and flavor transformation of supernova neutrinos. 

The numerical analysis is substantially simplified when using Eq. ( l49i) in Pjv approx- 
imation when N is at most as large as around ten. In P/v approximation the required 
computational power is proportional to numbers of equations:2(A^ + l)Ne where Nf, is the 
number of energy bins. In contrast, in the multi-angle simulation the number of equations 
in evolution is GNgNg where Ng is the number of angle bins and the factor 6 is the number 
of all types of neutrinos and anti-neutrinos. It is a huge amount of equations by noticing 
that NgjNg > 500 for L,^ = 10^^ erg/s [1]. More angle bins are required for larger neutrino 
luminosity in multi-angle simulation. The approach using moment equations reduces the 
computational power by about two orders of magnitude when is around ten. 

8 Numerical result 

In this section we show some numerical analysis on the non-linear flavor transformation 
of neutrinos. For simplicity we neglect matter effect in our analysis. 

The truncated set of Equations ( l50l) should be carefully treated in numerical analysis. 
Three terms in the right hand side of Eq. ( l50l) have different physics. The first term is 
caused by the diffusion of neutrinos and gives correction to the scaling law of moments. 
For example, Ql modifies the scaling of the total flux. The second term leaves Tr[py, 
Tr[(p'^)^] and Tr[{p'i^y] etc invariant. It is responsible for the unitary flavor evolution 
of moments. The third term does not change Tr[py but modifies Tr[(p'^)^], Tr[(p'^)^] 
etc. This term gives non-unitary evolution of neutrino flavor. The appearance of the 
non-unitary term in moment equations makes it hard to do numerical computations. 
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In the present analysis we simplify the numerical computation and search a self- 
consistent solution to the moment equations (l50l) . We assume that the evolution is 
dominated by the first and the second term in the right hand side of Eq. ( |50l) . We 
search for solution of the following form 

p', = + Apl (59) 

In this approximation Q^'^ and Ha,c are all expressed using 

Ql'" = Ql'\p'k'). HA,c = HA,c{Pk)- (60) 

satisfies 

^ =-ro'Ql{p^')-^[HA,p','] (61) 

For a small step we get 

p^r + At) = -—Qlipf) + e-^^-^>;°(r)e^^-^'- (62) 

^0 

Ap'^, is considered as a small correction and its initial condition is set to zero. It is sourced 

by p't- 

^ =-i[HA,Ap',]-i[Hc,Ql{pn- (63) 
Eq. ( l63i) can be re-written as 

^ = -^[H',^Ql{pl% (64) 

where 

Api = M\Ap',MA. pI' = M\p,'Ma, H'c = M\HcMa, (65) 
Ma = MA^Pk) is the evolution matrix given by Hamiltonian HA^Pk)'- 

j-MA = -iHaMa. (66) 

Initial condition Ma = 1 is set. For a small step Eq. (!64l) gives 

Ap, (r + Ar ) - Ap^ (r) = e'^^c Arg2^.Ht a. _ q2 ^q^) 

We expect Ap'^^ receive small contributions from p^ and its effect on the final result 
is negligible. Including Ap'^ in Q]f' and Ha^ introduces even smaller corrections. In 
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Figure 3: (color online) Fraction of ny^jijiy^ +^^1^2;)) versus radius r in two neutrino system. 
Two lines are for results with and without Ap^ correction separately. lArrig^^l = 3 x 10~^eV^, 
2 26*13 = 0.01. 
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numerical analysis the consistency of this approximation should be checked. That is, 
including Ap'^ should not give much modification to p^. 

This approximation is supported by the result of numerical computation. A result 
of numerical computation using this approximation is shown in Fig. [31 It is shown 
for two neutrino system of [yei'^^x) with inverted mass hierarchy and for emission angle 
distribution described by Eq. (ISSl) . We choose L^^ = Lp^ = L^^ = Lp^ = 3. x 10^^ erg/s. 
The initial energy spectrum of neutrino is given by the Fermi-Dirac distribution 

ME) = ,_f^. , (68) 

where x = E/T^, and F2 is the normalization factor. Parameters of four types of neutrinos 
are chosen as: T^, = 2.76 MeV, Tp^ = 4.01 MeV, T^, = Tp^ = 6.26 MeV. p^, = /ip, = 
f^fx = f^i^x = 3. In Fig. [3] one can see clearly the synchronized oscillation and the transition 
to bipolar oscillation. 

We show both result with Ap^ corrections and the result without Ap^ corrections. As 
can be seen these two results agree very well. Ap'^ corrections give negligible corrections 
to the neutrino evolution. This can be understood by noticing that in the small r region 
it is dominated by synchronized oscillation where all pf point to the same direction in 
flavor space and they commute with the effective Hamiltonian He- Hence effect of the 
third term in Eq. fISU]) does not change synchronized oscillation and is not important in 
small r region. In the large r region the magnitude of Ql is suppressed by rg/r^ and is 
again not important. This result shows that in this approximation we can neglect the 
third term in Eq. (150|) and use the following equation in our analysis 



dr ~ 



i[HA,Pk\ 



(69) 
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Figure 4: (color online) Fraction of ny^jijiy^ +^^1^2;)) versus radius r in two neutrino system. 
Two lines are for P4 and Pg approximations separately. Parameters of neutrinos are the same 
as in Fig. [3l 

Result in Fig. [3] is obtained using P4 approximation. In Fig. H] we compare the result 
of P4 approximation with the result of Pe approximation. One can see that results of 
these two approximations agree quite well. This result of numerical analysis shows that 
truncated moment equations indeed have very good convergence property. More detailed 
numerical analysis will be presented in other publications. 



9 Equations in two flavor system 

In two flavor system of neutrinos the density matrix can be expressed in terms of the unit 
matrix and the Pauli matrices. It's much easier to do analytic study in two flavor system. 
In this section we re-write the moment equations for two flavor system. 
For two flavors of neutrinos we write 



E'' P'AE) = \{P', +Pk-a), p[{E) = i(P^° + P^ ■ a). 



(70) 
(71) 



where E is the energy of neutrino and a = (cr^, cr^, a^) are Pauli matrices. 

Equations are obtained when using commutation relation of the Pauli matrices. Ne- 
glecting moments of > 1 and the third term in the right hand side of Eq. (150|) we get 
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in the mass base of neutrino 



dr 

dn 

dr 

dPi 
dr 

dP[ 



-2r^\-^Pi + {ueB + V2Gf{L + z'^Ni)] x P^, (72) 

-2r^\-^P[ + [~ujeB + V2Gf{L + r-^iVi)] x P^, (73) 

cjij^ + V2GfL + V2GF;z"^iVi) x A, (74) 

-ueB + V2G^L + V2Gfz-^Ni) x P{. (75) 



where oje = /S.m? /2E and 

Eqs. ( 1721) . ( 1731) . ( 17il) and ( 1751) form a closed set of equations. 

Eqs. (1721) . (173|) . (174|) and (175|) are equations we obtained for describing the pendulum 
motion of neutrinos in flavor space. Eqs. ([71]) and (175|) are the equations widely used 
in two flavor analysis in literature An improvement of Eqs ( l74l) and (|75l) is that the 
dependence on the radius r is explicitly derived. Correction to the scaling law of the 
Hamiltonian z~^N can be included when higher moments are included. This is another 
virtue using moment equations. 

We note that an important piece of this set of equations, Eqs. fl72l) and fl73l) . are 
missed in the literature. Actually what have been missed are more relevant to physical 
observables. Neutrino densities are described by Pq(Pq°) and Pq(Pq^). The dynamics of 
non-linear flavor transformation is controlled by Pi and P[ which are not real observables. 
Predictions of these two new equations on the flavor transformation should be explored. 

It should be reminded that Eqs. (1711) . fl7Sl) . (1721) and (175]) may not be reliably used in 
precise numerical analysis since moments of order k > 1 have all been neglected when 
deriving them. They can be used in qualitative analysis. We note that including the 
third term in the right hand side of Eq. (|50]) leads to extra term in pendulum equation 
of neutrinos. 



10 Conclusion 

In summary we derive a series of moment equations describing the transport and flavor 
transformation of neutrinos above neutrino sphere in core-collapse supernova. We examine 
the system of neutrinos in supernova using spherical coordinate and find a particular set of 
moments of the density matrix of neutrino. We expand the Liouville equation of neutrinos 
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using this series of moments and obtain the moment equations. We work with assumptions 
that the system is spherically symmetric and stationary in time scale around 10~^s. We 
arrive at a truncated set of moment equations (|50ll . After a further approximation we get 
(1691) . Rates of neutrino flavor transformation can be solved as functions of radius r using 
this set of truncated moment equations. These equations will be very helpful to future 
researches on supernova neutrinos [23]. 

We study geometric scaling properties of these moments in the expanding system of 
supernova and find that they have very good property of convergence. Higher moments 
Pk{k > 1) converge to zero much faster than lower moments po,i- Furthermore, we check 
initial conditions of moments and find that magnitudes of higher moments are naturally 
smaller than po,i on the neutrino sphere. Based on these observations we argue that the 
infinite series of moment equations can be safely truncated to be a small set of equations 
with pk{k < 10). Numerical analysis shows that the results in P/v approximation indeed 
converge for N ^ 10. In Fig. [3] and Fig. H] we show some numerical analysis and find 
patterns of synchronized oscillation and bipolar oscillation. 

The truncated set of moment equations have a number of good properties. First, 
the truncated set of equations tremendously simplify numerical analysis of the problem. 
Previous numerical works use many angle bins and energy bins to simulate the neutrino 
evolution. In multi-angle simulation evolutions of more than 500 trajectories have to be 
followed. Together with more than 500 energy bins and an extra factor 6 ( number of total 
types of neutrinos and anti- neutrinos) , more than a million equations have to be solved 
in this kind of simulation. It is ultra-complicated. The approach presented in this article 
uses a small set of truncated moments rather than many angle bins and is much simpler 
than the multi-angle simulation. Using a small set of truncated moment equations the 
required computational power in numerical analysis is reduced by two orders of magnitude 
compared to that in multi-angle simulation. 

Second, the moment equations systematically take into account the effect of the an- 
gular distribution of neutrino emission. As noted in Introduction, angular distribution of 
neutrino emission is essential for the effect of neutrino self-interaction to play an important 
in neutrino flavor transformation. The approach presented in the present article introduces 
moments of density matrix of neutrinos to describe the angular distribution of neutrinos. 
The equation of zeroth moment naturally includes the effect of higher moments on the 
evolution of neutrino intensity. The equations of higher moments systematically take into 
account the evolution of angular distribution of neutrino. In particular, modification to 
the scaling of the strength of the effective Hamiltonian can be included when including 
moments pk with k > 2. The dependence of this modification on the emission angle dis- 
tribution of neutrinos can be systematically studied. In contrast, previous researches use 
fixed scaling law (basically law) for the strength of the effective Hamiltonian. This set 
of moment equations provides a powerful tool towards a complete understanding to the 
effect of the emission angle distribution in flavor transformation of supernova neutrinos. 
More analysis on the effect of angular distribution in neutrino oscillation will be presented 
in other publications. 
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We also consider two flavor system of neutrinos. Using the truncated set of equations 
for zeroth and first moments we derive equations describing the pendulum motion of neu- 
trinos in flavor space. In addition to equations given in literature, we find new equations 
which are important and are more relevant to physical observable. 

There are still many problems. More detailed works should be done to check carefully 
the evolution of higher moments and their effect in convergence of moment equations. 
More numerical works should be done to analyze in detail the effect of higher moments on 
phenomena such as synchronized oscillation and bipolar oscillation, spectral split and adi- 
abaticity of neutrino evolution, etc. One important remaining problem is how to extend 
formulation to the case without the assumption of spherical symmetry. Answer to this 
problem will enable us to understand the impact of anisotropic disturbance of ordinary 
matter on neutrino flavor transformation. To answer this question we need to examine Eq. 
( IT2|) in more detail. Subtle problems in quantum evolution of moment equations should 
be investigated before we can answer this question. These problems will be explored in 
future works. 

Acknowledgement: I wish to thank Y. Z. Qian, G. Raffelt, A. Yu. Smirnov for discus- 
sions on neutrino flavor conversion in supernova. □ 



Appendix 

Neglecting the gravitational potential, the invariant distance in spherical coordinate = 
(r, 9, if) is given by 

ds^ = dt^ — gijdx^dx\ (77) 

where Qij is the metric 

gij = diag{l, r^, sin^ 9}. (78) 
The Christoffel symbol r*;^. is computed using the metric gij: 

We get 

Kb = -r, = -r sin' 9, T% = Vl = r~\ (80) 
= -sm9cos9, T^^ = V^^ = r-\ T^,^ = T^, = ctg9. (81) 
All other F's are zero. 



^Note added: After submission of this article I know from A. Mirizzi that an earher attempt to do 
moment expansion of the equation of supernova neutrinos is made in Ref. [25| . 
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Using Eq. (ITT!) we find 

dp' _ lr2(/)2 + r2sin2^(p^)2 
dt r p° 



(821 



It's easy to verify tliat 



— = ^ + sin^cos^ — ^, 83 

dt r p^ p" 

dp^ 2fp^ _ , „pV ...x 

-1- = — - 2 cigQ (84) 

dr r p^ p^ 



^Wf + r'(p^)' + sin^ ^(p^)^] = 0. (85) 



Tliis is just tlie statement tliat \p\ = \/'gijp^ = \/ {p'^Y + r'^{p^Y + sin^ 9{p'^y, the 
magnitude of the momentum, is a constant. It equals to p^ for neutrino with tiny masses. 

P^y y/Oeep' = ^p' and ^Jg^p"^ = r sin 9p'^ are the momenta projected to f, 9 and 
directions separately. Hence we can write 

p' = cos 9p p^, p' = - sin 9p cos fp p^, p^ = sin 9p sin ipp (86) 

r r sin 



Using Eq. ([82]) and ([83]) we get 



d cos 6*5, 1 , n ^ , , 

— -^ = - l-cos^^p , 87 
at r 



— ^ = — ct(76' sin 6'p sinyjp. (88) 
dt r 

The consequence of Eq. ( IHTI) is that cos 6'p approaches to 1. This is consistent with the 
physical picture that as neutrinos go out of the supernova their motion gets closer to the 
radial direction. 
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